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ABSTRACT 



v~j ■ High-resolution mid-infrared images have been obtained in N-band and Q- 

■ band for the proto-planetary nebula IRAS 16594—4656. A bright equatorial 

^ torus and a pair of bipolar lobes can clearly be seen in the infrared images. 

The torus appears thinner at the center than at the edges, suggesting that it is 
viewed nearly edge-on. The infrared lobes correspond to the brightest lobes of 
the reflection nebula seen in the Hubble Space Telescope (HST) optical image, 
but with no sign of the point-symmetric structure seen in the visible image. 
The lobe structure shows a close correspondence with a molecular hydrogen map 



^The paper is based on observations obtained at the Gemini Observatory. The Gemini Observatory is 
operated by the Association of Universities for Research in Astronomy, Inc., under a cooperative agree- 
ment with the NSF on behalf of the Gemini partnership: the National Science Foundation (United States), 
the Particle Physics and Astronomy Research Council (United Kingdom), the National Research Council 
(Canada), CONICYT (Chile), the Austrahan Research Council (Austraha), CNPq (Brazil) and CONICET 
(Argentina) . 
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obtained with HST^ suggesting that the dust emission in the lobes traces the 
distribution of the shocked gas. The shape of the bipolar lobes shows clearly 
that the fast outflow is still confined by the remnant circumstellar envelope of 
the progenitor asymptotic giant branch (AGB) star. However, the non-detection 
of the dust outside of the lobes suggests that the temperature of the dust in the 
AGB envelope is too low for it to be detected at 20 jim.. 

Subject headings: circumstellar matter: — infrared: stars — infrared: ISM: dust 
grains — ISM: planetary nebulae: general — stars: AGB and post-AGB 

1. INTRODUCTION 

Proto-planetary nebulae (PPNe) are the long-sought-after missing link between the end 
of the asymptotic giant branch (AGB) phase and the beginning of planetary nebula phase 
of stellar evolution. After the Infrared Astronomical Satellite {IRAS) mission, a number of 
objects were proposed as candidate PPNe based on their infrared colors and other spectral 

properties. These are typically stars of G to B spectral type with significant infrared excesses 
due to the remnant circumstellar dust shell ejected in the AGB phase. Of particular interest 
among these candidates arc a number of carbon-rich objects whose abundances show a strong 
enhancement of s-process elements, as expected from the dredge-up of material in thermal 
pulses during the AGB evolution (Kwok 1993; van Winckel 2003). For other candidates 
there is some possibility of confusion with massive supergiants, but those in this carbon-rich 
group are almost certainly bona-fide PPNe. 

One of these objects is IRAS 16594—4656. It is a bright mid-infrared source which has 
typical colors of a PPN (Volk & Kwok 1989). Optically it is associated with a southern 
emission-line star. It was found to be of spectral type B7 with V magnitude 14.6, subject 
to about 7.5 magnitudes of visual extinction (Van de Steene, van Hoof, & Wood 2000). It 
was not clear from the original IRAS spectral observations whether the object was oxygen- 
rich or carbon-rich, but subsequent Infrared Space Observatory {ISO) spectral observations 
showed that it has carbon-based dust features, including the 21 yum feature (Garcia-Lario et 
al. 1999). Optical images obtained with the Hubble Space Telescope showed that the star 
has a surrounding reflection nebula with a complex structure (Hrivnak, Kwok, & Su 1999). 
The relative faintness of the reflection nebula compared to the star, together with the optical 
morphology, led Hrivnak, Kwok, & Su (1999) to conclude that the nebula is intrinsically 
bipolar (or multi-polar) viewed at an intermediate angle to the bipolar axis. Both optical 
and near-infrared spectral observations show emission lines, which are thought to be shock- 
excited rather than radiatively excited by the star (Garcia-Lario et al. 1999; Van de Steene 
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& van Hoof 2003). No radio emission has been detected from IRAS 16594—4656 (Van de 
Steene & Pottasch 1993). 

The distance to this object is uncertain. The dust shell model of Hrivnak, Volk, & 
Kwok (2000) suggests a distance of 2.6 kpc if the total luminosity of the star is 10^ Lq. A 
slightly smaller estimate of (2.2±0.4) kpc is given by Van de Steene & van Hoof (2003) for 
the same assumed total luminosity. Most of the luminosity is being emitted in the infrared 
where extinction effects are smaller than in the optical, so these estimates are not strongly 
affected by the non-spherical morphology of the dust shell. 

The optical morphology of the nebula in IRAS 16594—4656 appears complex, with 
what appear to be pairs of symmetric structures at multiple position angles, for which it 
was named the "Water Lily Nebula" (Hrivnak, Kwok, & Su 1999). The optical images also 
showed several concentric arcs centered on the star (Hrivnak, Kwok, & Su 2001). However, 
near-infrared observations with the A^/CMO^" instrument on i/S'T showed a somewhat simpler 
morphology, although it was not clear whether this was simply an effect of reduced dynamic 
range compared to the optical observations (Su et al. 2001). An initial N-band observation of 
IRAS 16594-4656 with the TIMMI2 camera on the ESO 3.6m telescope failed to resolve the 
dust shell (Van de Steene, van Hoof, & Wood 2000), while more recent TIMMI2 observations 
did marginally resolve the structure (Garcia- Hernandez et al. 2004). 

In this paper we present new, higher sensitivity and higher angular resolution mid- 
infrared images of IRAS 16594—4656 which resolve the dust shell in thermal emission. We 
present the observations in section 2. We then derive a dust color temperature map in 
section 3, and compare the morphology as observed in the mid-infrared to that seen at other 
wavelengths in section 4. We give a brief discussion of the results in the final section. 

2. OBSERVATIONS 

The observations reported here were obtained with the T-ReCS instrument on Gemini 
South under program GS-2004A-Q-56. The sky and telescope background were removed 
by chopping and nodding during the observations. Images of IRAS 16594—4656 were ob- 
tained in three filters. In the N-band window: 360 second (total on-source exposure times) 
images with the "Si-5 11.66um" filter on 2004 March 11 and with the "Si-6 12.33um" fil- 
ter on 2004 May 8; and in Q-band window: a 600 second image with the "Qa 18.30um" 
filter on 2004 May 10. These filters will be referred to as the "Si5" , "Si6" , and "Qa" fil- 
ters, respectively. All these filters have widths A A ~ l//m. Information about these filters, 
including the filter profiles, can be found on the Gemini Observatory pubhc WWW pages 
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(see http://www.gemini.edu/sciops/instniments/miri/T-ReCSFiltersMtrnt). On each night, 
standard star observations for flux cahbration were carried out immediately after the obser- 
vation of the science target. These stars were HD 123139 on March 11, HD 169916 on May 
8, and HD 175775 on May 10. All of these stars are included among the mid-infrared spec- 
trophotometric standards of Cohen et al. (1999). Comparisons of the point-spread functions 
(PSFs) of these stars with observations of a Cen A/B or a CMa on various nights from 
December 2003 through May 2004 indicate that these stars are suitable as PSF references 
as well as spectral standards. 

The N-band observations were made under good conditions, as judged from the level 
of sky cancellation obtained while chopping. The Q-band observations were made under 
marginal conditions; however IRAS 16594—4656 is a very bright target at 18 /xm and it was 
detected with good signal-to-noise ratio despite the less than ideal conditions. The image 
quality was good as estimated from the standard stars. The full width at half maximum 
(FWHM) was 0'.'37 for the Si5 filter image, 0'.'39 for the Si6 image, and 0'.'60 for the Qa image. 
These correspond to Strehl values of about 0.6, fairly typical of better seeing conditions for 
the N-band filters but somewhat lower than usual for the Qa filter. 

The raw images on and off of the target were subtracted and then summed to produce 
raw images of IRAS 16594—4656 in the three filters. Flux calibration of the images was done 
in two different ways. First, the standard star observations were used to find the conversion 
from in-band counts to Jy using the assumed spectral energy distribution from Cohen et 
al. (1999) integrated over the filter profiles, and then these scaling factors were applied to 
the images of IRAS 16594—4656. The pixel by pixel brightnesses so obtained were then 
converted to Jy/square arc-second using the pixel size of T-ReCS. Second, the estimated 
filter flux densities were generated from the ISO spectrum of IRAS 16594—4656 in the same 
way as was used to get the expected flux densities in Jy for the standard stars. These values 
were also used to convert from counts to Jy/square arc-second in the images. It was found 
that for the Si5 fllter these two methods agreed within less than 1%. For the other fllters 
the agreement was poorer. Using the ISO spectrum for the Si6 filter estimate gave a value 
22% higher than that from the standard star. This discrepancy is too large to be due to 
an atmospheric extinction effect, judging from some estimates of the extinction coefficient 
in this filter made on other nights. It probably indicates that the sky conditions were not 
uniform in the directions to the standard star and to our target, as in other observations 
we have obtained with T-ReCS the inter-comparison with ISO spectra gives 2% agreement 
for filters in the N-band window. The Qa brightness calculated from the standard star 
came out 11% lower than that calculated from the ISO spectrum. This is probably within 
the uncertainties caused by the variable sky conditions. There is also a color effect due to 
differences in spectral shape between IRAS 16594—4656 and the standard stars, but since 
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the filters are relatively narrow this is a small correction and it was neglected. 

In what follows, we have chosen to use the ISO spectrum as the basis for creating surface 
brightness images, since this minimizes the effect of variable sky conditions. This assumes 
that the ISO spectrum gives the correct absolutely cahbrated total brightness and the T- 
ReCS observations give the correct relative brightness distribution, or equivalently that all 
the atmospheric effects are uniform over the small T-ReCS field of view. Figure 1 shows 
the three flux-calibrated images in Jy/square arc-second. The region shown in the figure 
is 7'.'2 X 7'.'2 and contains all the detected emission from IRAS 16594—4656. The total flux 
density for IRAS 16594-4656 was calculated to be 44.0, 56.9, and 177 Jy for the Si5, Si6, 
and Qa filters, respectively. 

In all three cases, the circumstellar shell of IRAS 16594—4656 appears as a bipolar 
nebula of dimension about 4'.'5 x 2'.'25, with a bright central region orientated roughly north- 
south and two lobes extending east and west. There is a clear difference in size between the 
east and west lobes. The east lobe is about 20% smaller than the west lobe both in width 
and maximum detected radius from the star in these images. The optical depth at these 
wavelengths, especially in Q-band, must be small; thus this size difference must be caused 
by either a physical size difference between the two lobes or distinctly different projection 
angles for the two lobes. 

The central bright region of the mid-infrared images appears to be some type of thin 
"equatorial" torus, perpendicular to the axis of the two lobes. The northern end of the this 
structure is brighter than the southern end in all the images, but the ratio is much closer to 
1:1 in the Qa image. This shows that the dust in the torus is cooler in the southern region 
than in the northern region. The sharp edge of the torus in the north is particularly striking, 
as the images are very bright there but there is a sudden edge beyond which no emission is 
observed along the line of the torus. 

3. IMAGE ANALYSIS 

3.1. Image Deconvolution 

Lucy deconvolution of the raw images was carried out using the standard star obser- 
vations as PSF templates. This was done using the stsdas.analysis.restore.lucy task in the 
STSCI reduction package under IRAF version 2.12a. A 61 by 61 pixel box was used to define 
the PSF. Pixels outside this PSF box but within a 181 by 181 pixel box centered on the star 
were used to derive the background level, which was subtracted from the stellar profile. 
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It was found that most of the improvement in the image resolution was obtained in the 
first few Lucy iterations, after which there was no significant change in the derived structure, 
so the deconvolution was stopped after 20 iterations. The deconvolved images were then re- 
smoothed with a gaussian of FWHM C'l, slightly larger than the original pixel size. The 
resulting images are sharpened by about a factor of 3 in PSF width compared to the original 
images. 

Figure 2 shows these sharpened images. The two N-band images have almost identical 
structure. The sharpened Q-band image has the same general morphology as the N-band 
images but both lobes arc seen to be smaller by about 0.4 arc-seconds than in the N-band 
images. However, we believe that this size difference is not real. Comparison of the raw 
images in the Qa and N-band filters shows that the emission region is just slightly larger for 
the Qa image than for the Si5 and Si6 filters. This indicates that the Lucy deconvolution 
for the lobes introduced a small artifact into the image. As the Strehl value for the standard 
star in the Qa filter was lower than usual, it is possible that the seeing changed between 
the target observation and the standard star observation. There was some indication of 
variable Q-band conditions during these observations. If the seeing did get worse for the 
standard star observation, that would explain the decrease in lobe size for the deconvolved 
image compared to the raw image, although the magnitude of the decrease is larger than 
one would expect based upon the FWHM value for the standard star Qa filter image. 

The structure of the brightest regions is similar for all three filters. The bright part 
suggests a torus of some sort, and it appears to be thinner at the center than at the edges. 
This suggests that the torus is seen nearly edge on. If it were oriented at some intermediate 
angle to the plane of the sky, as was earlier asserted based upon the visible images, then 
one would expect the torus to appear as a small ellipse in these mid-infrared images. This 
is clearly not seen. However, if we are indeed viewing the torus edge-on then the visible and 
mid-infrared images suggest it to be quite asymmetric with position. 

In each panel of Figure 2, the estimated star position is marked with a small black dot. 
This position was found by cross-comparing the T-ReCS Si5 image with an HST NICMOS 
image taken in the narrow-band filter centered on the H2 2.12 //m line (Hrivnak, Kelly, & 
Su 2004). The estimated stellar position is very close to the geometrical center of the bright 
"bar" in the dust emission. It is located in a region of relatively low brightness in the three 
filters, especially in the Si5 and Si6 filters. This probably means that in the optical images 
we are seeing the star through some type of hole in the toroid, where little dust is present. 



- 7- 



3.2. Dust Color Temperature Maps 



From the surface brightness images in two fihers it is possible to construct a ratio map, 
and if one assumes that the surface brightness is due to optically-thin thermal emission from 
dust grains of a known type, then the surface brightness ratio can be converted into color 
temperatures between the two wavelengths. Denoting the surface brightness in Jy/square 
arc-second by S{i'), the color temperature, Tc, is defined by solving 



for each brightness ratio value in the image. The optical depth values Tj, for the two filters 
are assumed to be proportional to the absorption cross-sections, so we can replace the optical 
depths with the actual Qabs values for this calculation. The scattering component of the 
dust extinction is expected to be small at these long wavelengths. The ratio of r^, values is 
a constant as long as the dust grain properties are uniform in the dust shell, so then the 
equation gives a one-to-one transformation between surface brightness ratio and Tc. 

We have used the Si5 and Qa images to produce a brightness ratio map, after convolving 
the Si5 image with a gaussian profile of FWHM of 0'.'195 to match the angular resolution 
of the Qa image. Regions of the two images which had "low" brightness values, taken to 
be less than 0.1% of the respective peak values, were masked out of the ratio image. The 
ratio image was then transformed to dust color temperature, assuming that the dust grains 
are 0.1 /^m amorphous carbon grains with the opacity function for AC type 2 grains from 
Rouleau & Martin (1991). These grain properties were the basis of the spectral model for 
IRAS 16594—4656 presented by Hrivnak, Volk, & Kwok (2000). For another assumed dust 
grain size or type the dust color temperature values would change, but the relative variations 
over the image should remain the same. While the values do not indicate the physical 
temperatures of the dust grains, since they are some type of average along the hne of sight 
for each pixel, they do indicate the global dust temperature variations in the circumstellar 
shell as long as the dust grains are not drastically different than assumed for the calculation. 
The Tc map is shown in the lower right panel of Figure 1. The Tc values are confined to a 
relatively narrow range. The bright part of the dust shell has a range of Tc from about 140 
K to 160 K. It was found that the Tc map was much the same whether or not the Si5 image 
was convolved to the resolution of the Qa image. 

The Tc map shows that the region of highest dust color temperature is in the north part 
of the central bright region, while the color temperature is much lower directly to the south 
on the other side of the stellar position. There is also a region of high dust color temperature 
in the southern wall of the west lobe. 



T^-,By{Tc, Vl) 
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There is a cap at the end of the east lobe that is seen in the color temperature map, 
which is also apparent in the Si5 image. This cap does not appear to be associated with a 
higher dust color temperature than elsewhere in lobe, which suggests that the dust optical 
depth is higher here than for other positions in the lobe. It is possible that the east lobe is 
smaller than the west lobe because its expansion is impeded by external material, and that 
the cap represents a boundary between the lobe and the external medium. 

4. DISCUSSION and CONCLUSIONS 

Figure 3 shows the deconvolved Si5 image with overlaid contour plots in the optical 
I-band (0.8 /xm; Su, Hrivnak, & Kwok 2001) and ncar-infrarcd H2 filter images (2.12 ^m; 
Hrivnak, Kelly, & Su 2004). This allows direct comparison of the morphology of the nebula 
in these different wavelengths. The H2 image was matched to the Si5 filter image since 
they were immediately seen to have very similar morphologies. This was used to obtain 
a reasonably accurate estimate of the central star position in the T-ReCS image. That 
position was then used to match the optical image to the mid-infrared image. The Si5 image 
is plotted in absolute surface brightness units, Jy/square arc-second. 

As shown in the upper right panel of Figure 3, there is a close match of features in the 
lobes between the NICMOS H2 filter image and the T-ReCS image. The bright lobe edges 
are regions of strong H2 emission. The cap in the east lobe is clearly visible in this panel, and 
it also corresponds to a region of strong H2 emission. It is more difficult to determine if any 
of the mid-infrared bright waist structure is also detected in the H2 image, because the star 
is saturated in the II2 image, but it does not look as if anything but the edges of the lobes 
are detected in the II2 image. Since the II2 emission is mainly shock excited (Van de Stccnc 
& van Hoof 2003), this raises the possibility that the dust emission from these regions are 
partially shock-excited. Another possibility that could explain why the dust emission region 
so closely matches the shock is that the dust grains may be much smaller downstream from 
the shock than they are before the shock, and that as a result the small grains are transiently 
heated to relatively high temperatures. 

Comparison of the T-ReCS images with the iKT optical images (Hrivnak, Kwok, & Su 
1999) shows that the two lobes seen in the T-ReCS images correspond closely to the central 
brightest part of the reflection nebula. The total size of the optical reflection nebula is 12.3" 
by 8.8", which is much larger than the size of the N-band or Q-band images. The optical 
image overlay in the lower left panel of Figure 3 uses logarithmically spaced contours ranging 
from 0.0085% of the stellar peak brightness up to the peak brightness, with each contour 
at a level 2.7 times the previous contour. The optical reflection nebula is at most about 
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0.75% of the stellar peak brightness. The lowest contours show what has been suggested to 
be point-symmetric morphology, with three pairs of oppositely directed features at position 
angles of about 40°, 57°, and 87° east of north as can be measured from the I-band image 
presented in Hrivnak, Kwok, & Su (1999). The optical lobes have been suggested to be 
caused by a rapidly preccssing, columnatcd high-speed wind from the star (Garcia-Lario 
et al. 1999). The T-ReCS image shows just the two lobes at position angle 75°. While 
there seems to be some very faint mid-infrared emission detected in the Si5 filter on a size 
scale roughly three times that of the main mid-infrared lobes, over-plotting this with the 
optical image does not show any correspondences with the faint extended optical structure. 
In particular, none of this emission is seen just outside the bright torus either to the north or 
to the south. Kinematic observations of the individual optical lobes is needed to determine 
if they are distinct structures or not. If these point symmetric features are not distinct 
kinematic structures, then perhaps there are just two lobes but the optical appearance is 
due to structure in the walls of these lobes which make them look more complex. 

The optical emission is due to reflection by dust, and so the larger size of the optical 
images compared to the mid-infrared images shows that the dust shell is much larger than 
is obvious from the N-band or Q-band images. The dust outside the central bright region 
delineated by the shocked H2 emission is clearly much colder than that inside the shocks. 
At the ends of mid-infrared torus, in particular, there must be a large discontinuity in tem- 
perature and optical depth so that there is a large change in mid-infrared surface brightness 
but a much smaller change in the optical brightness of the scattered hght. 

These mid-infrared images suggest that the torus is perpendicular to the plane of the 
sky. This is consistent with recent observations of the bipolar lobes, which indicate that they 
are nearly in the plane of the sky. Unpublished high resolution long-slit spectra in the near- 
infrared (Hrivnak et al. 2006) have been obtained with the Phoenix spectrograph on Gemini 
South. Analysis of these spectra, which map the molecular hydrogen line at 2.12 fj,m for 
cuts at different position angles through the star and the lobes, shows similar velocities for 
the two lobes, indicating that they are oriented very close to the plane of the sky. This was 
also concluded by Ueta, Murakawa, & Meixner (2005) based upon near-infrared polarimetry 
along with dust shell modeling of the spectral energy distribution of the object. Since the 
evidence indicates that the nebular axis is oriented very close to the plane of the sky, the 
east lobe must actually be smaller than the west lobe. 

The T-ReCS images show that the central bright region of the circumstellar dust shell 
is quite different than that of other well-known bipolar nebulae IRAS 17150—3224, IRAS 
17441—2411, Roberts 22, and Hen 3-401, in all of which the dust emission is strongly peaked 
at the stellar position. For IRAS 16594—4656 the mid-infrared brightness is at a minimum 
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near the stellar position and is much higher to the north and the south along the central 
waist. Possibly the torus is of much lower optical depth in this object than in the others, or it 
is highly asymmetric with a low optical depth along our line of sight. The latter suggestion is 
consistent with the visibility of the central star. Certainly for IRAS 17150—3224 and IRAS 
17441—2411 the spectral energy distribution implies a relatively high optical depth along 
our line of sight to the star, and the star is not seen in visible light. The spectral type of 
the star in IRAS 16594—4656 is much earlier than that in the latter two objects, so it may 
simply be more evolved and the torus may have had more time to disperse. 

Unlike the hour-glass or open lobes observed in most bipolar planetary nebulae (e.g. 
NGC 6302), the bipolar lobes of IRAS 16594—4656 as shown in Figure 2 are closed and 
resemble the lobes seen in the PPN IRAS 17106—3046 and the young planetary nebula 
Hen 2-320. The morphology of the lobes clearly shows that the lobes are confined by the 
circumstellar medium, and the fast coUimated outflow which creates the bipolar lobes has 
not yet broken through the stellar wind of the AGB progenitor. The interaction between 
the fast and slow winds is clearly delineated by dust distribution in the lobes. We also note 
that there is "bulge" at the tip of the western lobe, which suggests that the high- velocity 
flow is on the verge of breaking out. In contrast to the cap at the tip of the eastern lobe 
(which represents a pile-up at the wind interface), the western lobe may represent a slightly 
more advanced stage of the breakout, and therefore explains the difference in sizes between 
the two lobes. In a few hundred years, we expect that both lobes will open up into butterfly 
morphology. We are therefore witnessing a critical phase of morphological transformation of 
PNs. 

We have successfully detected the bipolar lobes and a central bright waist in mid- 
infrared images of IRAS 16594—4656. The bright waist suggests that we are seeing a central 
torus nearly edge-on. While this is consistent with published polarization and unpublished 
kinematic results, it differs from earlier published interpretations of the visible image that 

concluded that the lobes are seen at an intermediate orientation. This emphasizes the 
need for multi-wavelength observations to confidently understand the structure of proto- 
planetary nebulae. This result may well be applicable to our understanding of other bipolar 
phenomenon such as YSOs and AGNs. 

B. J.H acknowledges support by the National Science Foundation under Grant No. 0407087. 
This work was supported in part by grants to S.K from the Natural Sciences and Engineering 
Research Council of Canada. 
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Fig. 1. — The flux-calibrated T-ReCS false-color images of IRAS 16594—4656. The four 
panels give respectively the Si5 image (upper left), the Si6 image (upper right), the Qa image 
(lower left), and the color temperature map derived from the Si5 to Qa surface brightness 
ratio (lower right). The T-ReCS images are given as surface brightnesses in units of Jy/square 
arc-second. The color temperature values are given in K. For each panel the bar gives the 
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color mapping. The image sections are all 7'.'17 square. All images have north up and east 
to the left. 
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Fig. 2. — The deconvolved T-ReCS false-color images of IRAS 16594—4656. The panels 
correspond to those in Figure 1. In these images the estimated position of the star is marked 
by a small dot. 
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Fig. 3. — A comparison of the morphology of IRAS 16594—4656 observed at optical, near- 
infrared, and mid-infrared wavelengths. The upper left panel shows the deconvolved Si5 
image as in Figure 2. The field of view is 14'.'34 by 10'.'75. In the lower left panel the same 
figure is shown with logarithmically spaced contours from the HST I band image of (Su et 
al. 2001) overlaid. The upper right panel shows the Si5 filter image again, with square-root 
scaling to show the low level emission better; this panel is magnified by a factor of 2 to 
show the inner quadrant of the left panels. Overlaid on the upper right T-ReCS image are 
linear contours from an H2 image obtained with the HST NICMOS narrow-band 2.1212 //m 
emission-line filter (Hrivnak, Kelly, & Su 2004). 



